Analysis MA and correlation
# ##### Create empty dataframe with slopes and confidence interval for all pairwise
Estimates_pairwise <- data.frame(Variables = rep(rep(c("correlation","cor_CI_inf","cor_CI_sup",
"slope", "intercept"),2),3),
Generation = rep(c(rep(7, 5), rep(29, 5)),3),
Pairwise = rep(c("Cherry_Cranberry",
"Cranberry_Strawberry",
"Strawberry_Cherry"), each=10),
Estimates = NA)
#########################
######################### G7
#########################
#######
####### Cherry Cranberry
#######
## Compute weighted correlation
weightedcor <- sjstats:::weighted_correlation(TEMP_dataG7_CheCran, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)
## Compute Major Axis regression
model <- lmodel2::lmodel2(logchange_allop ~ logchange_symp,
range.y = "interval",range.x = "interval",
data = TEMP_dataG7_CheCran[rep(1:nrow(TEMP_dataG7_CheCran), TEMP_dataG7_CheCran$N_sumsympallop),], nperm=0)
## The estimation of the confidence interval with cor.test is wrong
cor.test(x=TEMP_dataG7_CheCran[rep(1:nrow(TEMP_dataG7_CheCran), TEMP_dataG7_CheCran$N_sumsympallop), "logchange_symp"], y=TEMP_dataG7_CheCran[rep(1:nrow(TEMP_dataG7_CheCran), TEMP_dataG7_CheCran$N_sumsympallop), "logchange_allop"])
##
## Pearson's product-moment correlation
##
## data: TEMP_dataG7_CheCran[rep(1:nrow(TEMP_dataG7_CheCran), TEMP_dataG7_CheCran$N_sumsympallop), and TEMP_dataG7_CheCran[rep(1:nrow(TEMP_dataG7_CheCran), TEMP_dataG7_CheCran$N_sumsympallop), "logchange_symp"] and "logchange_allop"]
## t = 6.3642, df = 136, p-value = 2.79e-09
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.3390707 0.5982491
## sample estimates:
## cor
## 0.4790334
# Fill in data
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"] <- weightedcor$estimate
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- weightedcor$ci[1]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- weightedcor$ci[2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="intercept"] <- model$regression.results[2, 2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="slope"] <- model$regression.results[2, 3]
#######
####### Cranberry Strawberry
#######
## Compute weighted correlation
weightedcor <- sjstats:::weighted_correlation(TEMP_dataG7_CranStraw, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)
## Compute Major Axis regression
model <- lmodel2::lmodel2(logchange_allop ~ logchange_symp,
range.y = "interval",range.x = "interval",
data = TEMP_dataG7_CranStraw[rep(1:nrow(TEMP_dataG7_CranStraw), TEMP_dataG7_CranStraw$N_sumsympallop),], nperm=0)
# Fill in data
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="correlation"] <- weightedcor$estimate
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- weightedcor$ci[1]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- weightedcor$ci[2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="intercept"] <- model$regression.results[2, 2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="slope"] <- model$regression.results[2, 3]
#######
####### Strawberry Cherry
#######
## Compute weighted correlation
weightedcor <- sjstats:::weighted_correlation(TEMP_dataG7_StrawChe, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)
## Compute Major Axis regression
model <- lmodel2::lmodel2(logchange_allop ~ logchange_symp,
range.y = "interval",range.x = "interval",
data = TEMP_dataG7_StrawChe[rep(1:nrow(TEMP_dataG7_StrawChe), TEMP_dataG7_StrawChe$N_sumsympallop),], nperm=0)
# Fill in data
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="correlation"] <- weightedcor$estimate
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- weightedcor$ci[1]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- weightedcor$ci[2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="intercept"] <- model$regression.results[2, 2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="slope"] <- model$regression.results[2, 3]
#########################
######################### G29
#########################
#######
####### Cherry Cranberry
#######
## Compute unweighted correlation
unweightedcor <- cor.test(x=TEMP_dataG29_CheCran$logchange_symp, y=TEMP_dataG29_CheCran$logchange_allop)
## Compute weighted correlation
weightedcor <- sjstats:::weighted_correlation(TEMP_dataG29_CheCran, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)
## Compute Major Axis regression
model <- lmodel2::lmodel2(logchange_allop ~ logchange_symp,
range.y = "interval",range.x = "interval",
data = TEMP_dataG29_CheCran[rep(1:nrow(TEMP_dataG29_CheCran), TEMP_dataG29_CheCran$N_sumsympallop),], nperm=0)
# Fill in data
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"] <- unweightedcor$estimate
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- unweightedcor$conf.int[1]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- unweightedcor$conf.int[2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="intercept"] <- model$regression.results[2, 2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="slope"] <- model$regression.results[2, 3]
#######
####### Cranberry Strawberry
#######
## Compute unweighted correlation
unweightedcor <- cor.test(x=TEMP_dataG29_CranStraw$logchange_symp, y=TEMP_dataG29_CranStraw$logchange_allop)
## Compute weighted correlation
weightedcor <- sjstats:::weighted_correlation(TEMP_dataG29_CranStraw, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)
## Compute Major Axis regression
model <- lmodel2::lmodel2(logchange_allop ~ logchange_symp,
range.y = "interval",range.x = "interval",
data = TEMP_dataG29_CranStraw[rep(1:nrow(TEMP_dataG29_CranStraw), TEMP_dataG29_CranStraw$N_sumsympallop),], nperm=0)
# Fill in data
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="correlation"] <- unweightedcor$estimate
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- unweightedcor$conf.int[1]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- unweightedcor$conf.int[2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="intercept"] <- model$regression.results[2, 2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="slope"] <- model$regression.results[2, 3]
#######
####### Strawberry Cherry
#######
## Compute unweighted correlation
unweightedcor <- cor.test(x=TEMP_dataG29_StrawChe$logchange_symp, y=TEMP_dataG29_StrawChe$logchange_allop)
## Compute weighted correlation
weightedcor <- sjstats:::weighted_correlation(TEMP_dataG29_StrawChe, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)
## Compute Major Axis regression
model <- lmodel2::lmodel2(logchange_allop ~ logchange_symp,
range.y = "interval",range.x = "interval",
data = TEMP_dataG29_StrawChe[rep(1:nrow(TEMP_dataG29_StrawChe), TEMP_dataG29_StrawChe$N_sumsympallop),], nperm=0)
# Fill in data
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="correlation"] <- unweightedcor$estimate
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- unweightedcor$conf.int[1]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- unweightedcor$conf.int[2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="intercept"] <- model$regression.results[2, 2]
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="slope"] <- model$regression.results[2, 3]
Estimates_pairwise
## Variables Generation Pairwise Estimates
## 1 correlation 7 Cherry_Cranberry 0.479033353
## 2 cor_CI_inf 7 Cherry_Cranberry 0.001431407
## 3 cor_CI_sup 7 Cherry_Cranberry 0.956635299
## 4 slope 7 Cherry_Cranberry 4.160858972
## 5 intercept 7 Cherry_Cranberry -0.356412153
## 6 correlation 29 Cherry_Cranberry -0.567557085
## 7 cor_CI_inf 29 Cherry_Cranberry -0.908773300
## 8 cor_CI_sup 29 Cherry_Cranberry 0.228504360
## 9 slope 29 Cherry_Cranberry -1.427127889
## 10 intercept 29 Cherry_Cranberry 0.268405728
## 11 correlation 7 Cranberry_Strawberry 0.535360085
## 12 cor_CI_inf 7 Cranberry_Strawberry 0.105487696
## 13 cor_CI_sup 7 Cranberry_Strawberry 0.965232474
## 14 slope 7 Cranberry_Strawberry 1.070761415
## 15 intercept 7 Cranberry_Strawberry -0.049525449
## 16 correlation 29 Cranberry_Strawberry 0.593494172
## 17 cor_CI_inf 29 Cranberry_Strawberry -0.191100790
## 18 cor_CI_sup 29 Cranberry_Strawberry 0.915350061
## 19 slope 29 Cranberry_Strawberry 1.613908561
## 20 intercept 29 Cranberry_Strawberry -0.582936853
## 21 correlation 7 Strawberry_Cherry 0.616674051
## 22 cor_CI_inf 7 Strawberry_Cherry 0.188354905
## 23 cor_CI_sup 7 Strawberry_Cherry 1.044993196
## 24 slope 7 Strawberry_Cherry 0.829196290
## 25 intercept 7 Strawberry_Cherry 0.177499018
## 26 correlation 29 Strawberry_Cherry 0.676000585
## 27 cor_CI_inf 29 Strawberry_Cherry -0.300322260
## 28 cor_CI_sup 29 Strawberry_Cherry 0.960575095
## 29 slope 29 Strawberry_Cherry 2.537639762
## 30 intercept 29 Strawberry_Cherry -0.764708220
Compute confidence interval of the difference
computeCIcordifG29G7(pair="Cherry_Cranberry")
## dif lower higher
## -1.0465904 -1.9126980 -0.1804829
computeCIcordifG29G7(pair="Cranberry_Strawberry")
## dif lower higher
## 0.05813409 -0.78991100 0.90617917
computeCIcordifG29G7(pair="Strawberry_Cherry")
## dif lower higher
## 0.05932653 -0.95762427 1.07627734
# ##### Create empty dataframe with slopes and confidence interval for all pairwise
Estimates_pairwise <- data.frame(Variables = rep(rep(c("correlation","cor_CI_inf","cor_CI_sup",
"slope", "intercept"),2),3),
Generation = rep(c(rep(7, 5), rep(29, 5)),3),
Pairwise = rep(c("Cherry_Cranberry",
"Cranberry_Strawberry",
"Strawberry_Cherry"), each=10),
Estimates = NA ,
EstimatesWeighed = NA ,
EstimatesUnweighed = NA )
# Calculate sd in case of NA value
mean_sd_G7 <- sqrt(mean(data_logchange$sd_logchange[data_logchange$Generation=="7"]^2, na.rm = TRUE))
mean_sd_G29 <- sqrt(mean(data_logchange$sd_logchange[data_logchange$Generation=="29"]^2, na.rm = TRUE))
## Check whether the squared standard errors are of the same magnitude across fruit media (answer:yes)
tapply(data_logchange$sd_logchange^2, list(data_logchange$Treatment, data_logchange$Generation), mean, na.rm = TRUE)
## Impute directly the missing standard errors
data_logchange$se_logchange <- data_logchange$sd_logchange
data_logchange$se_logchange[data_logchange$Generation=="7"&is.na(data_logchange$se_logchange)] <- rep(mean_sd_G7, sum(data_logchange$Generation=="7"&is.na(data_logchange$se_logchange)))
## Nbr sample
nb_sim = 5000
#########################
######################### G7
#########################
#######
####### Cherry Cranberry
#######
# MA axis regression with bootstrap to estimate slope, intercept and
sim <- data.frame(t(sapply(1:nb_sim, bootstrap_lmodel2, data_pairwise = TEMP_dataG7_CheCran, sd_NA = mean_sd_G7)))
# Fill out data with quantile of simulations
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"] <- mean(sim$correlation.cor, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- quantile(sim$correlation.cor, probs = 0.025, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- quantile(sim$correlation.cor, probs = 0.975, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="intercept"] <- mean(sim$intercept, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="slope"] <- mean(sim$slope, na.rm = TRUE)
## Add other correlations
Estimates_pairwise$EstimatesWeighed[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"] <- weighted_correlation(TEMP_dataG7_CheCran, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)$estimate
Estimates_pairwise$EstimatesWeighed[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- weighted_correlation(TEMP_dataG7_CheCran, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)$ci[1]
Estimates_pairwise$EstimatesWeighed[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- weighted_correlation(TEMP_dataG7_CheCran, x=logchange_symp, y=logchange_allop, weights=N_sumsympallop, ci.lvl = 0.95)$ci[2]
Estimates_pairwise$EstimatesUnweighed[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"] <- cor.test(x=TEMP_dataG7_CheCran$logchange_symp, y=TEMP_dataG7_CheCran$logchange_allop)$estimate
Estimates_pairwise$EstimatesUnweighed[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- cor.test(x=TEMP_dataG7_CheCran$logchange_symp, y=TEMP_dataG7_CheCran$logchange_allop)$conf.int[1]
Estimates_pairwise$EstimatesUnweighed[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- cor.test(x=TEMP_dataG7_CheCran$logchange_symp, y=TEMP_dataG7_CheCran$logchange_allop)$conf.int[2]
#Save sim G7 to calculate proportion later:
sim_CheCranG7<-sim
#######
####### Cranberry Strawberry
#######
# MA axis regression with bootstrap to estimate slope, intercept and
sim <- data.frame(t(sapply(1:nb_sim, bootstrap_lmodel2, data_pairwise = TEMP_dataG7_CranStraw, sd_NA = mean_sd_G7)))
# Fill out data with quantile of simulations
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="correlation"] <- quantile(sim$cor, probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- quantile(sim$cor, probs = 0.025, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- quantile(sim$cor, probs = 0.975, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="intercept"] <- quantile(sim$intercept,probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="slope"] <- quantile(sim$slope,probs = 0.5, na.rm = TRUE)
## Save sim G7 to calculate correlation difference later:
sim_CranStrawG7<-sim
#######
####### Strawberry Cherry
#######
# MA axis regression with bootstrap to estimate slope, intercept and
sim <- data.frame(t(sapply(1:nb_sim, bootstrap_lmodel2, data_pairwise = TEMP_dataG7_StrawChe, sd_NA = mean_sd_G7)))
# Fill out data with quantile of simulations
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="correlation"] <- quantile(sim$cor, probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- quantile(sim$cor, probs = 0.025, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- quantile(sim$cor, probs = 0.975, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="intercept"] <- quantile(sim$intercept,probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==7&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="slope"] <- quantile(sim$slope,probs = 0.5, na.rm = TRUE)
## Save sim G7 to calculate correlation difference later:
sim_StrawCheG7<-sim
#########################
######################### G29
#########################
#######
####### Cherry Cranberry
#######
# MA axis regression with bootstrap to estimate slope, intercept and
sim <- data.frame(t(sapply(1:nb_sim, bootstrap_lmodel2, data_pairwise = TEMP_dataG29_CheCran, sd_NA = mean_sd_G29)))
# Fill out data with quantile of simulations
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"] <- quantile(sim$cor, probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- quantile(sim$cor, probs = 0.025, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- quantile(sim$cor, probs = 0.975, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="intercept"] <- quantile(sim$intercept, probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="slope"] <- quantile(sim$slope, probs = 0.5, na.rm = TRUE)
## Save sim G29 to calculate correlation difference later:
sim_CheCranG29<-sim
#######
####### Cranberry Strawberry
#######
# MA axis regression with bootstrap to estimate slope, intercept and
sim <- data.frame(t(sapply(1:nb_sim, bootstrap_lmodel2, data_pairwise = TEMP_dataG29_CranStraw, sd_NA = mean_sd_G29)))
# Fill out data with quantile of simulations
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="correlation"] <- quantile(sim$cor, probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- quantile(sim$cor, probs = 0.025, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- quantile(sim$cor, probs = 0.975, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="intercept"] <- quantile(sim$intercept,probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cranberry_Strawberry"&
Estimates_pairwise$Variables=="slope"] <- quantile(sim$slope,probs = 0.5, na.rm = TRUE)
## Save sim G29 to calculate correlation difference later:
sim_CranStrawG29<-sim
#######
####### Strawberry Cherry
#######
rm(sim)
# MA axis regression with bootstrap to estimate slope, intercept and
sim <- data.frame(t(sapply(1:nb_sim, bootstrap_lmodel2, data_pairwise = TEMP_dataG29_StrawChe, sd_NA = mean_sd_G29)))
# Fill out data with quantile of simulations
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="correlation"] <- quantile(sim$cor, probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_inf"] <- quantile(sim$cor, probs = 0.025, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="cor_CI_sup"] <- quantile(sim$cor, probs = 0.975, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="intercept"] <- quantile(sim$intercept,probs = 0.5, na.rm = TRUE)
Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Strawberry_Cherry"&
Estimates_pairwise$Variables=="slope"] <- quantile(sim$slope,probs = 0.5, na.rm = TRUE)
Estimates_pairwise
## Save sim G29 to calculate correlation difference later:
sim_StrawCheG29<-sim
## Compute the confidence interval of the difference between correlations in G7 and G29
difcorCheCran <- sim_CheCranG29$correlation.cor-sim_CheCranG7$correlation.cor
quantile(difcorCheCran, probs = c(0.025, 0.975), na.rm=TRUE)
difcorCranStraw <- sim_CranStrawG29$correlation.cor-sim_CranStrawG7$correlation.cor
quantile(difcorCranStraw, probs = c(0.025, 0.975), na.rm=TRUE)
difcorStrawChe <- sim_StrawCheG29$correlation.cor-sim_StrawCheG7$correlation.cor
quantile(difcorStrawChe, probs = c(0.025, 0.975), na.rm=TRUE)
#Cherry Cranbery
#Compute proportion of simulations correlations G7 are lower than correlation estimated on G29
cor_G29_CheCran<-Estimates_pairwise$Estimates[Estimates_pairwise$Generation==29&
Estimates_pairwise$Pairwise=="Cherry_Cranberry"&
Estimates_pairwise$Variables=="correlation"]
#Estimates of correlation during final phenotyping for Cherry Cranberry
p_val_decreaseCheCran<-length(sim_CheCranG7$correlation.cor[sim_CheCranG7$correlation.cor<cor_G29_CheCran])/nb_sim
p_val_decreaseCheCran
Plot
ymin = -50
ymax = 50
## Function equation R
eq_r <- function(gen = 7, pair = "Cherry_Cranberry") {
if(gen == 7 &pair == "Cherry_Cranberry"){
as.character(
as.expression(
substitute(~~italic(r)^2~"="~r2~"["~inf~","~sup~"]",
list(r2 = format(Estimates_pairwise$Estimates[Estimates_pairwise$Generation == gen &
Estimates_pairwise$Pairwise == pair &
Estimates_pairwise$Variables == "correlation"], digits = 2),
inf = format(Estimates_pairwise$Estimates[Estimates_pairwise$Generation == gen &
Estimates_pairwise$Pairwise == pair &
Estimates_pairwise$Variables == "cor_CI_inf"], digits = 1, nsmall=2),
sup = format(Estimates_pairwise$Estimates[Estimates_pairwise$Generation == gen &
Estimates_pairwise$Pairwise == pair &
Estimates_pairwise$Variables == "cor_CI_sup"], digits = 2)))
)
)
}else{
as.character(
as.expression(
substitute(~~italic(r)^2~"="~r2~"["~inf~","~sup~"]",
list(r2 = format(Estimates_pairwise$Estimates[Estimates_pairwise$Generation == gen &
Estimates_pairwise$Pairwise == pair &
Estimates_pairwise$Variables == "correlation"], digits = 2),
inf = format(Estimates_pairwise$Estimates[Estimates_pairwise$Generation == gen &
Estimates_pairwise$Pairwise == pair &
Estimates_pairwise$Variables == "cor_CI_inf"], digits = 2, nsmall=2),
sup = format(Estimates_pairwise$Estimates[Estimates_pairwise$Generation == gen &
Estimates_pairwise$Pairwise == pair &
Estimates_pairwise$Variables == "cor_CI_sup"], digits = 2)))
)
)
}
}
##################################################
################## G7
# Limits
ymin_CheCranG7=min(min(TEMP_dataG7_CheCran$logchange_allop-1.96*TEMP_dataG7_CheCran$sd_allop, na.rm= TRUE),
min(TEMP_dataG7_CheCran$logchange_symp-1.96*TEMP_dataG7_CheCran$sd_symp, na.rm= TRUE))
ymax_CheCranG7=max(max(TEMP_dataG7_CheCran$logchange_allop + 1.96*TEMP_dataG7_CheCran$sd_allop, na.rm= TRUE),
max(TEMP_dataG7_CheCran$logchange_symp + 1.96*TEMP_dataG7_CheCran$sd_symp, na.rm= TRUE))
lim_text<-ymin_CheCranG7+0.99*(ymax_CheCranG7-ymin_CheCranG7)
# Plot
CheCran_G7 <- ggplot(data = TEMP_dataG7_CheCran) +
geom_abline(intercept = Estimates_pairwise$Estimates[Estimates_pairwise$Generation == 7 &
Estimates_pairwise$Pairwise == "Cherry_Cranberry" &
Estimates_pairwise$Variables == "intercept"],
slope = Estimates_pairwise$Estimates[Estimates_pairwise$Generation == 7 &
Estimates_pairwise$Pairwise == "Cherry_Cranberry" &
Estimates_pairwise$Variables == "slope"],
colour = "black", size = 0.75) +
geom_text(x = -0.825, y = lim_text, label = eq_r(gen = 7, pair = "Cherry_Cranberry"), parse = TRUE, color="black", size = 3.5) +
geom_errorbar(aes(x = logchange_symp, ymin = logchange_allop-(1.96*sd_allop),
ymax = logchange_allop + (1.96*sd_allop),
color = Symp, linetype = Line_type),
width= 0.02, size = 0.3, alpha = 0.8) +
geom_errorbarh(aes(y = logchange_allop, xmin = logchange_symp-1.96*sd_symp, xmax = logchange_symp + 1.96*sd_symp,
color = Symp, linetype = Line_type),
height = 0.02, size = 0.3, alpha = 0.8) +
geom_point(aes(x = logchange_symp, y = logchange_allop, color = Symp,fill = Symp, shape = Allop),
size = 3, fill = "white", stroke =1.2) +
xlab("Fitness change in\nselective environment") +
ylab("Fitness change in\nalternative environment") +
ggtitle("Cherry vs. Cranberry") +
theme_LO_sober +
labs(shape = "Test fruit", color = "Evolution fruit") +
scale_shape_manual(values = c(21, 22)) +
scale_color_manual(values = c("#BC3C6D", "#FDB424")) +
coord_cartesian(ylim = c(ymin_CheCranG7, ymax_CheCranG7),
xlim = c(ymin_CheCranG7, ymax_CheCranG7))
CheCran_G7

# Limits
ymin_CranStrawG7=min(min(TEMP_dataG7_CranStraw$logchange_allop-1.96*TEMP_dataG7_CranStraw$sd_allop, na.rm= TRUE),
min(TEMP_dataG7_CranStraw$logchange_symp-1.96*TEMP_dataG7_CranStraw$sd_symp, na.rm= TRUE))
ymax_CranStrawG7=max(max(TEMP_dataG7_CranStraw$logchange_allop + 1.96*TEMP_dataG7_CranStraw$sd_allop, na.rm= TRUE),
max(TEMP_dataG7_CranStraw$logchange_symp + 1.96*TEMP_dataG7_CranStraw$sd_symp, na.rm= TRUE))
lim_text<-ymin_CranStrawG7+0.99*(ymax_CranStrawG7-ymin_CranStrawG7)
CranStraw_G7 <- ggplot(data = TEMP_dataG7_CranStraw) +
geom_abline(intercept = Estimates_pairwise$Estimates[Estimates_pairwise$Generation == 7 &
Estimates_pairwise$Pairwise == "Cranberry_Strawberry" &
Estimates_pairwise$Variables == "intercept"],
slope = Estimates_pairwise$Estimates[Estimates_pairwise$Generation == 7 &
Estimates_pairwise$Pairwise == "Cranberry_Strawberry" &
Estimates_pairwise$Variables == "slope"],
colour = "black", size = 0.75) +
geom_errorbar(aes(x =logchange_symp, ymin = logchange_allop-(1.96*sd_allop),
ymax = logchange_allop + (1.96*sd_allop),
color = Symp, linetype = Line_type),
width= 0.02, size = 0.3, alpha = 0.8) +
geom_errorbarh(aes(y = logchange_allop, xmin = logchange_symp-1.96*sd_symp, xmax = logchange_symp + 1.96*sd_symp,
color = Symp, linetype = Line_type),
height = 0.02, size = 0.3, alpha = 0.8) +
geom_point(aes(x =logchange_symp, y = logchange_allop, color = Symp,fill = Symp, shape = Allop),
size =3, fill = "white", stroke =1.2) +
geom_text(x = -0.45, y = lim_text, label = eq_r(gen = 7, pair = "Cranberry_Strawberry"), parse = TRUE, color="black", size = 3.5) +
xlab("Fitness change in\nselective environment") +
ylab("Fitness change in\nalternative environment") +
ggtitle("Cranberry vs. Strawberry") +
coord_cartesian(ylim = c(ymin_CranStrawG7, ymax_CranStrawG7),
xlim = c(ymin_CranStrawG7, ymax_CranStrawG7)) +
labs(shape = "Test fruit", color = "Selection fruit") +
scale_shape_manual(values = c(22, 24)) +
scale_color_manual(values = c("#FDB424", "#3FAA96")) +
theme_LO_sober
CranStraw_G7

# Limits
ymin_StrawCheG7=min(min(TEMP_dataG7_StrawChe$logchange_allop-1.96*TEMP_dataG7_StrawChe$sd_allop, na.rm= TRUE),
min(TEMP_dataG7_StrawChe$logchange_symp-1.96*TEMP_dataG7_StrawChe$sd_symp, na.rm= TRUE))
ymax_StrawCheG7=max(max(TEMP_dataG7_StrawChe$logchange_allop + 1.96*TEMP_dataG7_StrawChe$sd_allop, na.rm= TRUE),
max(TEMP_dataG7_StrawChe$logchange_symp + 1.96*TEMP_dataG7_StrawChe$sd_symp, na.rm= TRUE))
lim_text<-ymin_StrawCheG7+0.99*(ymax_StrawCheG7-ymin_StrawCheG7)
StrawChe_G7 <- ggplot(data = TEMP_dataG7_StrawChe) +
geom_abline(intercept = Estimates_pairwise$Estimates[Estimates_pairwise$Generation == 7 &
Estimates_pairwise$Pairwise == "Strawberry_Cherry" &
Estimates_pairwise$Variables == "intercept"],
slope = Estimates_pairwise$Estimates[Estimates_pairwise$Generation == 7 &
Estimates_pairwise$Pairwise == "Strawberry_Cherry" &
Estimates_pairwise$Variables == "slope"],
colour = "black", size = 0.75) +
geom_errorbar(aes(x =logchange_symp, ymin = logchange_allop-(1.96*sd_allop), ymax = logchange_allop + (1.96*sd_allop),
color = Symp, linetype = Line_type),
width= 0.02, size = 0.3, alpha = 0.8) +
geom_errorbarh(aes(y = logchange_allop, xmin = logchange_symp-1.96*sd_symp, xmax = logchange_symp + 1.96*sd_symp,
color = Symp, linetype = Line_type),
height = 0.02, size = 0.3, alpha = 0.8) +
geom_point(aes(x = logchange_symp, y = logchange_allop, color = Symp,fill = Symp, shape = Allop),
size =3, fill = "white", stroke =1.2) +
geom_text(x = -0.83, y = lim_text, label = eq_r(gen = 7, pair = "Strawberry_Cherry"), parse = TRUE, color="black", size = 3.5) +
xlab("Fitness change in\nselective environment") +
ylab("Fitness change in\nalternative environment") +
ggtitle("Strawberry vs. Cherry") +
labs(shape = "Test fruit", color = "Selection fruit") +
scale_shape_manual(values = c(21, 24)) +
scale_color_manual(values = c("#BC3C6D", "#3FAA96")) +
coord_cartesian(ylim = c(ymin_StrawCheG7, ymax_StrawCheG7),
xlim = c(ymin_StrawCheG7, ymax_StrawCheG7)) +
theme_LO_sober
StrawChe_G7

#####################################################################################
################################## G29 ##################################
#####################################################################################
TEMP_dataG29_CheCran
## Line Treatment Fruit_s Generation N_allop logchange sd_logchange SA
## 15 CEA Cranberry Cherry 29 30 -0.21774071 0.13281291 0
## 17 CEB Cranberry Cherry 29 30 -0.59200204 0.09634940 0
## 21 CEC Cranberry Cherry 29 30 -0.20670829 0.12591289 0
## 37 CRA Cherry Cranberry 29 30 0.32537794 0.10025669 0
## 41 CRB Cherry Cranberry 29 30 0.16929862 0.12699030 0
## 45 CRC Cherry Cranberry 29 30 0.08474123 0.12334598 0
## 47 CRD Cherry Cranberry 29 30 0.54736475 0.09718018 0
## 49 CRE Cherry Cranberry 29 30 -0.10496044 0.11556266 0
## Symp Allop Line_type logchange_allop sd_allop logchange_symp
## 15 Cherry Cranberry dashed -0.21774071 0.13281291 0.14767000
## 17 Cherry Cranberry dashed -0.59200204 0.09634940 0.25038029
## 21 Cherry Cranberry dashed -0.20670829 0.12591289 0.35854062
## 37 Cranberry Cherry dashed 0.32537794 0.10025669 -0.14885350
## 41 Cranberry Cherry dashed 0.16929862 0.12699030 0.29719689
## 45 Cranberry Cherry dashed 0.08474123 0.12334598 0.02703717
## 47 Cranberry Cherry dashed 0.54736475 0.09718018 -0.15627744
## 49 Cranberry Cherry dashed -0.10496044 0.11556266 0.72513486
## sd_symp N_symp N_sumsympallop
## 15 0.1574284 30 60
## 17 0.1009864 30 60
## 21 0.1169300 30 60
## 37 0.1193552 30 60
## 41 0.1531197 30 60
## 45 0.1156381 30 60
## 47 0.1156695 30 60
## 49 0.1275430 30 60
TEMP_dataG29_CranStraw
## Line Treatment Fruit_s Generation N_allop logchange sd_logchange SA
## 39 CRA Strawberry Cranberry 29 30 -0.29731422 0.1098141 0
## 42 CRB Strawberry Cranberry 29 30 -0.24539416 0.1387043 0
## 43 CRC Strawberry Cranberry 29 30 -0.28950778 0.2012393 0
## 46 CRD Strawberry Cranberry 29 30 -1.02555272 0.1902082 0
## 51 CRE Strawberry Cranberry 29 30 -0.16810249 0.1398057 0
## 67 FRA Cranberry Strawberry 29 30 -0.25158762 0.1712755 0
## 72 FRB Cranberry Strawberry 29 30 0.04303751 0.1441119 0
## 73 FRC Cranberry Strawberry 29 30 0.52109773 0.1254820 0
## Symp Allop Line_type logchange_allop sd_allop logchange_symp
## 39 Cranberry Strawberry dashed -0.29731422 0.1098141 -0.14885350
## 42 Cranberry Strawberry dashed -0.24539416 0.1387043 0.29719689
## 43 Cranberry Strawberry dashed -0.28950778 0.2012393 0.02703717
## 46 Cranberry Strawberry dashed -1.02555272 0.1902082 -0.15627744
## 51 Cranberry Strawberry dashed -0.16810249 0.1398057 0.72513486
## 67 Strawberry Cranberry dashed -0.25158762 0.1712755 0.15709032
## 72 Strawberry Cranberry dashed 0.04303751 0.1441119 0.56023100
## 73 Strawberry Cranberry dashed 0.52109773 0.1254820 0.36640738
## sd_symp N_symp N_sumsympallop
## 39 0.1193552 30 60
## 42 0.1531197 30 60
## 43 0.1156381 30 60
## 46 0.1156695 30 60
## 51 0.1275430 30 60
## 67 0.1411560 30 60
## 72 0.1284283 30 60
## 73 0.1201170 30 60
TEMP_dataG29_StrawChe
## Line Treatment Fruit_s Generation N_allop logchange sd_logchange SA
## 14 CEA Strawberry Cherry 29 30 0.01971359 0.1475176 0
## 16 CEB Strawberry Cherry 29 30 -0.15850680 0.0983504 0
## 19 CEC Strawberry Cherry 29 30 0.21131464 0.1237821 0
## 68 FRA Cherry Strawberry 29 30 -0.50599095 0.1378652 0
## 71 FRB Cherry Strawberry 29 30 0.19377147 0.1094819 0
## 75 FRC Cherry Strawberry 29 30 0.32151693 0.1140574 0
## Symp Allop Line_type logchange_allop sd_allop logchange_symp
## 14 Cherry Strawberry dashed 0.01971359 0.1475176 0.1476700
## 16 Cherry Strawberry dashed -0.15850680 0.0983504 0.2503803
## 19 Cherry Strawberry dashed 0.21131464 0.1237821 0.3585406
## 68 Strawberry Cherry dashed -0.50599095 0.1378652 0.1570903
## 71 Strawberry Cherry dashed 0.19377147 0.1094819 0.5602310
## 75 Strawberry Cherry dashed 0.32151693 0.1140574 0.3664074
## sd_symp N_symp N_sumsympallop
## 14 0.1574284 30 60
## 16 0.1009864 30 60
## 19 0.1169300 30 60
## 68 0.1411560 30 60
## 71 0.1284283 30 60
## 75 0.1201170 30 60
# Limits
ymin_CheCranG29=min(min(TEMP_dataG29_CheCran$logchange_allop-1.96*TEMP_dataG29_CheCran$sd_allop, na.rm= TRUE),
min(TEMP_dataG29_CheCran$logchange_symp-1.96*TEMP_dataG29_CheCran$sd_symp, na.rm= TRUE))
ymax_CheCranG29=max(max(TEMP_dataG29_CheCran$logchange_allop + 1.96*TEMP_dataG29_CheCran$sd_allop, na.rm= TRUE),
max(TEMP_dataG29_CheCran$logchange_symp + 1.96*TEMP_dataG29_CheCran$sd_symp, na.rm= TRUE))
lim_text<-ymin_CheCranG29+0.99*(ymax_CheCranG29-ymin_CheCranG29)
CheCran_G29 <- ggplot(data = TEMP_dataG29_CheCran) +
geom_abline(intercept = Estimates_pairwise$Estimates[Estimates_pairwise$Generation ==29 &
Estimates_pairwise$Pairwise == "Cherry_Cranberry" &
Estimates_pairwise$Variables == "intercept"],
slope = Estimates_pairwise$Estimates[Estimates_pairwise$Generation ==29 &
Estimates_pairwise$Pairwise == "Cherry_Cranberry" &
Estimates_pairwise$Variables == "slope"],
colour = "black", size = 0.75) +
geom_errorbar(aes(x =logchange_symp, ymin = logchange_allop-(1.96*sd_allop),
ymax = logchange_allop + (1.96*sd_allop),
color = Symp),
width= 0.02, size = 0.3, alpha = 0.8) +
geom_errorbarh(aes(y = logchange_allop, xmin = logchange_symp-1.96*sd_symp, xmax = logchange_symp + 1.96*sd_symp,
color = Symp),
height = 0.02, size = 0.3, alpha = 0.8) +
geom_point(aes(x =logchange_symp, y = logchange_allop, color = Symp,fill = Symp, shape = Allop),
size =3, fill = "white", stroke =1.2) +
geom_text(x = 0.33, y = lim_text, label = eq_r(gen = 29, pair = "Cherry_Cranberry"), parse = TRUE, color="black", size = 3.5) +
xlab("Fitness change in\nselective environment") +
ylab("Fitness change in\nalternative environment") +
ggtitle("Cherry vs. Cranberry") +
labs(shape = "Test fruit", color = "Evolution fruit") +
scale_shape_manual(values = c(16, 15)) +
scale_color_manual(values = c("#BC3C6D", "#FDB424")) +
coord_cartesian(ylim = c(ymin_CheCranG29, ymax_CheCranG29),
xlim = c(ymin_CheCranG29, ymax_CheCranG29)) +
theme_LO_sober
CheCran_G29

# Limits
ymin_CranStrawG29=min(min(TEMP_dataG29_CranStraw$logchange_allop-1.96*TEMP_dataG29_CranStraw$sd_allop, na.rm= TRUE),
min(TEMP_dataG29_CranStraw$logchange_symp-1.96*TEMP_dataG29_CranStraw$sd_symp, na.rm= TRUE))
ymax_CranStrawG29=max(max(TEMP_dataG29_CranStraw$logchange_allop + 1.96*TEMP_dataG29_CranStraw$sd_allop, na.rm= TRUE),
max(TEMP_dataG29_CranStraw$logchange_symp + 1.96*TEMP_dataG29_CranStraw$sd_symp, na.rm= TRUE))
lim_text<-ymin_CranStrawG29+0.99*(ymax_CranStrawG29-ymin_CranStrawG29)
#Plot
CranStraw_G29 <- ggplot(data = TEMP_dataG29_CranStraw) +
geom_abline(intercept = Estimates_pairwise$Estimates[Estimates_pairwise$Generation ==29 &
Estimates_pairwise$Pairwise == "Cranberry_Strawberry" &
Estimates_pairwise$Variables == "intercept"],
slope = Estimates_pairwise$Estimates[Estimates_pairwise$Generation ==29 &
Estimates_pairwise$Pairwise == "Cranberry_Strawberry" &
Estimates_pairwise$Variables == "slope"],
colour = "black", size = 0.75) +
geom_errorbar(aes(x =logchange_symp, ymin = logchange_allop-(1.96*sd_allop),
ymax = logchange_allop + (1.96*sd_allop),
color = Symp),
width= 0.02, size = 0.3, alpha = 0.8) +
geom_errorbarh(aes(y = logchange_allop, xmin = logchange_symp-1.96*sd_symp, xmax = logchange_symp + 1.96*sd_symp,
color = Symp),
height = 0.02, size = 0.3, alpha = 0.8) +
geom_point(aes(x =logchange_symp, y = logchange_allop, color = Symp,fill = Symp, shape = Allop),
size =3, fill = "white", stroke =1.2) +
geom_text(x = -0.45, y = lim_text, label = eq_r(gen = 29, pair = "Cranberry_Strawberry"), parse = TRUE, color="black", size = 3.5) +
xlab("Fitness change in\nselective environment") +
ylab("Fitness change in\nalternative environment") +
ggtitle("Cranberry vs. Strawberry") +
coord_cartesian(ylim = c(ymin_CranStrawG29, ymax_CranStrawG29),
xlim = c(ymin_CranStrawG29, ymax_CranStrawG29)) +
labs(shape = "Test fruit", color = "Selection fruit") +
scale_shape_manual(values = c(15, 17)) +
scale_color_manual(values = c("#FDB424", "#3FAA96")) +
theme_LO_sober
CranStraw_G29

# Limits
ymin_StrawCheG29=min(min(TEMP_dataG29_StrawChe$logchange_allop-1.96*TEMP_dataG29_StrawChe$sd_allop, na.rm= TRUE),
min(TEMP_dataG29_StrawChe$logchange_symp-1.96*TEMP_dataG29_StrawChe$sd_symp, na.rm= TRUE))
ymax_StrawCheG29=max(max(TEMP_dataG29_StrawChe$logchange_allop + 1.96*TEMP_dataG29_StrawChe$sd_allop, na.rm= TRUE),
max(TEMP_dataG29_StrawChe$logchange_symp + 1.96*TEMP_dataG29_StrawChe$sd_symp, na.rm= TRUE))
lim_text<-ymin_StrawCheG29+0.99*(ymax_StrawCheG29-ymin_StrawCheG29)
StrawChe_G29 <- ggplot(data = TEMP_dataG29_StrawChe) +
geom_abline(intercept = Estimates_pairwise$Estimates[Estimates_pairwise$Generation ==29 &
Estimates_pairwise$Pairwise == "Strawberry_Cherry" &
Estimates_pairwise$Variables == "intercept"],
slope = Estimates_pairwise$Estimates[Estimates_pairwise$Generation ==29 &
Estimates_pairwise$Pairwise == "Strawberry_Cherry" &
Estimates_pairwise$Variables == "slope"],
colour = "black", size = 0.75) +
geom_errorbar(aes(x =logchange_symp, ymin = logchange_allop-(1.96*sd_allop),
ymax = logchange_allop + (1.96*sd_allop),
color = Symp),
width= 0.02, size = 0.2, alpha =1) +
geom_errorbarh(aes(y = logchange_allop, xmin = logchange_symp-1.96*sd_symp, xmax = logchange_symp + 1.96*sd_symp,
color = Symp),
height = 0.02, size = 0.2, alpha =1) +
geom_point(aes(x =logchange_symp, y = logchange_allop, color = Symp,fill = Symp, shape = Allop),
size =3, fill = "white", stroke =1.2) +
geom_text(x = -0.2, y = lim_text, label = eq_r(gen = 29, pair = "Strawberry_Cherry"), parse = TRUE, color="black", size = 3.5) +
coord_cartesian(ylim = c(ymin_StrawCheG29, ymax_StrawCheG29),
xlim = c(ymin_StrawCheG29, ymax_StrawCheG29)) +
xlab("Fitness change in\nselective environment") +
ylab("Fitness change in\nalternative environment") +
ggtitle("Strawberry vs. Cherry") +
labs(shape = "Test fruit", color = "Selection fruit") +
scale_shape_manual(values = c(16, 17)) +
scale_color_manual(values = c("#BC3C6D", "#3FAA96")) +
theme_LO_sober
StrawChe_G29

legend_tradevide <- ggplot(data = data_logchange[data_logchange$Generation == "29",],
aes(x =logchange, y = logchange, color = Symp,fill = Symp, shape = Allop)) +
geom_point(size =2.5, fill = "white") +
labs(shape = "Test fruit", color = "Selection fruit") +
scale_shape_manual(values = c(16, 15, 17)) +
scale_color_manual(values = c("#BC3C6D", "#FDB424", "#3FAA96")) +
theme_LO_sober
legend_trade <- lemon::g_legend(legend_tradevide)
Slopeestimates_logchange_pairwise_errorbar <- cowplot::ggdraw() +
cowplot::draw_plot(CheCran_G7 + theme(legend.position = "none"),
x = 0.01, y = 0.5, width = 0.22, height = 0.45) +
cowplot::draw_plot(CranStraw_G7 + theme(legend.position = "none"),
x = 0.31, y = 0.5, width = 0.22, height = 0.45) +
cowplot::draw_plot(StrawChe_G7 + theme(legend.position = "none"),
x = 0.61, y = 0.5, width = 0.22, height = 0.45) +
cowplot::draw_plot(legend_trade, x = 0.85, y = 0.5, width = 0.1, height = 0.1) +
cowplot::draw_plot(CheCran_G29 + theme(legend.position = "none"),
x = 0.01, y = 0, width = 0.22, height = 0.45) +
cowplot::draw_plot(CranStraw_G29 + theme(legend.position = "none"),
x = 0.31, y = 0, width = 0.22, height = 0.45) +
cowplot::draw_plot(StrawChe_G29 + theme(legend.position = "none"),
x = 0.61, y = 0, width = 0.22, height = 0.45) +
cowplot::draw_plot_label(c("Intermediate phenotyping step", "A", "B", "C", " ",
"Final phenotyping step", "D", "E", "F", " "),
x = c(0.30, 0.01, 0.30, 0.61, 0.92, 0.250, 0.01, 0.30, 0.61, 0.92),
y = c(0.99, 0.95, 0.95, 0.95, 0.95, 0.49, 0.45, 0.45, 0.45, 0.45),
hjust = c(-0.25, -0.25, -0.25, -0.25, -0.25, -0.75, -0.75, -0.75, -0.75, -0.75),
vjust = c(1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5),
size = 16)
Slopeestimates_logchange_pairwise_errorbar

#
# #
#
cowplot::save_plot(file =here::here("figures", "FIG_Correlation.pdf"),
Slopeestimates_logchange_pairwise_errorbar,
base_height = 18/cm(1), base_width = 34/cm(1), dpi = 610)